Draft version December 3, 2010 

Preprint typeset using 1^1^^ style emulatcapj v. 03/07/07 



A MULTI-EPOCH STUDY OF THE RADIO CONTINUUM EMISSION OF ORION SOURCE I: CONSTRAINTS 
ON THE DISK EVOLUTION OF A MASSIVE YSO AND THE DYNAMICAL HISTORY OF ORION BN/KL 

C. GoDDi, E. M. L. Humphreys 

European Southern Observatory, Karl-Schwarzschild-Strasse 2, D-85748 Garching bei Mjinchen and 
Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138 



o 

(N 
o 

Q 

o 

> 

0^ 
l> 
m 



X 



L. J. Greenhill 

Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138 

AND 

C. J. Chandler 

National Radio Astronomy Observatory, P.O. Box O, Socorro, NM 87801 

AND 

L. D. Matthews 

MIT Haystack Observatory, Westford, MA 01886 
Draft version December 3, 2010 

ABSTRACT 

We present new A7 mm continuum observations of Orion BN/KL with the Very Large Array (VLA). 
We resolve the emission from the Young Stehar Objects (YSO) radio Source I and BN at several epochs. 
Radio Source I is highly elongated northwest-southeast, and remarkably stable in flux density, position 
angle, and overall morphology over nearly a decade. This favors the extended emission component 
arising from an ionized edge-on disk rather than an outwardly propagating jet. 

We have measured the proper motions of Source I and BN for the first time at 43 GHz. We 
confirm that both sources are moving at high speed (12 and 26 km s~^ , respectively) approximately 
in opposite directions, as previously inferred from measurements at lower frequencies. We discuss 
dynamical scenarios that can explain the large motions of both BN and Source I and the presence of 
disks around both. Our new measurements support the hypothesis that a close (~ 50 AU) dynamical 
interaction occurred around 500 years ago between Source I and BN as proposed by Gomez et al. From 
the dynamics of encounter we argue that Source I today is likely to be a binary with a total mass on 
the order of 20 M©, and that it probably existed as a softer binary before the close encounter. This 
enables preservation of the original accretion disk, though truncated to its present radius of ^ 50 AU. 

N-body numerical simulations show that the dynamical interaction between a binary of 20 M© total 
mass (Source I) and a single star of 10 M© mass (BN) may lead to the ejection of both and binary 
hardening. The gravitational energy released in the process would be large enough to power the 
wide-angle, high- velocity flow traced by H2 and CO emission in the BN/KL nebula. Assuming the 
proposed dynamical history is correct, the smaller mass for Source I recently estimated from SiO 
maser dynamics (> 7 M©) by Matthews et al., suggests that non-gravitational forces (e.g. magnetic) 
must play an important role in the circumstellar gas dynamics. 

Subject headings: ISM: individual objects: Orion BN/KL - stars: formation - (stars:) binaries (includ- 
ing multiple): close - methods: N-body simulations 



1. INTRODUCTION 

The Orion BN/KL complex, at a d istance of 418 ±6 pc 
(|Menten et al.ll2007t iKim et al.l l2008'). contains the near- 
est region of ongoing high-mass star formation. A dense 
protostellar cluster lies within the region containing three 
radio sources that are believed to be massive young stel- 
lar objects (YSOs): the highly embedded radio Source I, 
(iGreenhill et al.' 2004a: Rcid ct al. 2007: Matt hews et alJ 
I2OIOD : the BN object, which is the brightest source in the 
regio n in the mid-infrared (IR) at 12.4 /xm (jCezari et al.l 
Il998| l: and Source n, a relativ ely evolved YSO w ith a 
disk observed in the MIR (Gr eenhill et al.l l2004bD and 
a jet observed in the rad io at 8.4 GHz ()Menten &: ReidI 
[IQM [G6mez et al.l[2008h . 

Despite intensive investigations at radio and IR wave- 
lengths, the primary heating source(s) for the Orion KL 



region [L ^ 10^ L©) is (are) still not known. An- 
other long-standing puzzle is the geometry of outflow 
and the identification of driving sources. There are 
two large-scale outflows in the region. A powerful 
(3 X lO'^'^ ergs), high-velocity (30-200 km s'^ ), wide- 
angle (^ 1 rad) outflow extends northwest-southeast 
(NW-SE) over 0.3 pc. This so-call ed "high-velocity" 
outflow is traced in CO emission dCh ernin fc Wright! 
119961 : IZapata et aDl2009t ) and in 2.12 ^m H2 shocked 
emission originating in finger-likc structures that end 
in bow shocks ( Allen fc Burton .1993) . A second, "low- 
velocity" outflow (~ 18 km s~^ ) is identified by a clus- 
ter of bright w = SiO and H2O masers, concen- 
trated within a few arcsec around S ource I and elon- 
gated northeast-southwest (NE-SW; iGenzel fc Stutzkl 
Il989.: Gree nhill ct al. 1998 , and in prep.). 

Source I has been proposed as a possible driver of 
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both the high- velocity NW-SE outflow (e.g., Wright et al. 
1995: lGreenhill et alJfToM iBallv fc Zinneckei]l2005ll and 



the low-velocity NE-SW outflow (jBeuther et al.l 120051: 
iGreenhiU et"alTl2004al and in prep.). Confusion arises 
because the radio continuum emission from Source I 
shows an elongated structure, which has been inter- 
preted as both an i onize d jet along a NW-SE direction 
i|Ballv fc Zinneckeij[200l lT^ l2008a^ and as a n ionized 
disk with a NE-SW spin axis (jReid et al.ll2007[ ). 

Based on a multi-epoch observational campaign of sev- 
eral SiO maser transitions using the Very Large Ar- 
ray (VLA) and t he Ve ry Long Baseline Array (VLB A), 
iMatthews et al.l (|2010l ) and Greenhill et al. (in prep.) 
provide convincing evidence that Source I is associated 
with a disk/outflow system with a NE-SW axis. In par- 
ticular, based on a VL B A monitoring oiv = 1, 2 J = 1—0 
SiO maser transitions, IMatthews et al.l ()2010D presented 
a movie of bulk gas flow tracing a compact disk and 
the base of a protostellar outflow at radii < 100 AU 
from Source I. In addition, Greenhill et al. measured 
proper motions of w = SiO masers, which trace the 
bipolar outflow expanding with a characteristic velocity 
~ 18 km s~^ and extending to radii of 100-1000 AU from 
Source I along a NE-SW axis, parallel to the axis of the 
disk/outflow system traced by w = 1,2 masers at radii 
< 100 AU. 

The origin and nature of the wide-angle NW-SE ori- 
ented outflow, traced in s hocked H2 and CO emi ssion, 
is still a matter of debate. iRodriguez et al.l ()2005| ) pro- 
posed that Source I and BN had a close encounter about 
500 yrs ago, which resulted in the ejection of interacting 
sources and the formation of a tight binary (Source I). 
Based on a proper motion s tudy of radio sources at 
8.4 GHz, iGomez et all ()2008l ) proposed that Source I, 
BN, and source n participated in the same dynamical in- 
teraction 500 yrs ago and that all three sources are mov- 
ing away from the putative center of interaction. The 
energy liberated in the process may provide in princi- 
ple sufficient energy to p ower the fast NW-SE outflow 
(iBallv fc Zinneckeill2005t IZapata et ani2009( l. It is not 
clear, however, what effect a close passage and conse- 
quent formation of a tight binary would have on a well- 
organized ac cretion /outflow struct ur e such as observed 
in Source I (jMatthews et al.l[20Tl . iTa^ ([2003 . 12008H ) 
proposed an alternatively scenario where a close pas- 
sage between Source I and BN (a runaway star from the 
Trapezium) would trigger a tidally-enhanced accretion 
and subsequent outburst of outflow activity, resulting in 
the powerful high-velocity outflow. 

In this paper, we present new multi-epoch, high- 
angular resolution observations of the radio continuum 
emission in Orion BN/KL at 7 mm (43 GHz) from the 
VLA. The main goals of the new observations were to 
investigate the nature of the radio continuum in Source I 
and reconstruct the dynamical history of BN/KL. In par- 
ticular, we investigate changes in morphology, size, and 
flux density as a function of time (over a decade) and 
frequency (43 vs 8.4 GHz) to obtain insights into the na- 
ture of the radio sources, mainly to test the ionized disk 
and jet hypotheses proposed for Source I. In addition, 
we measured absolute proper motions of radio sources 
based on accurate absolute astrometry, with the aim of 
constraining the dynamical history of the BN/KL region. 
In order to quantify probabilities of different dynami- 



cal scenarios, we present also new N-body simulations 
of decaying protostellar clusters with a varying num- 
ber of objects . Previous N-body simulations for BN/KL 
()G6mez et al.l l2008) assumed a five-member cluster hav- 
ing large masses (in the range 8-20 Mq), which resulted 
in the formation of binaries with total mass of 36 Mq. 
However, there is no evidence of such massive objects 
in the BN/KL region, based on present data. Our new 
simulations assume more plausible masses of individual 
objects as well as investigate a larger number of possible 
scenarios. 

The current paper is structured as follows. The obser- 
vational setup and data calibration procedures are de- 
scribed in § 2. § [3] reports the results of the multi-epoch 
study. In § |4l we discuss the morphological evolution 
of the radio continuum from Source I and its interpre- 
tation in terms of an ionized disk. In § |S] we suggest 
that Source I and BN had a past close passage, based on 
proper motion measurements. In § |6l we discuss dynam- 
ical scenarios that can explain the origin of the proper 
motions measured for Source I and BN. § 7 and 8 discuss 
problems related to the mass of Source I and a possible 
origin for the fast, wide-angle H2 outflow, respectively. 
Finally, conclusions are drawn in § 9. 

2. OBSERVATIONS AND DATA REDUCTION 

Observations of the A7 mm continuum emission in 
Orion BN/KL were made using the Very Large Array 
(VLA) of the National Radio Astronomy Observatory 
(NRAO)^. We present data obtained in four distinct 
epochs spanning over 8 years (see Table 1). All data 
were obtained while the VLA was in the A-configuration, 
yielding a resolution of approximately 0"05^. 

In order to image the continuum emission from radio 
sources in BN/KL and to provide a strong signal as a 
phase-reference to calibrate the phase and amplitude of 
the (weak) continuum signal, we employed a dual-band 
continuum setup with a narrow band (6.25 MHz in pro- 
grams AC952 and AC817; 3.25 MHz in AG622; and 
1.56 MHz in AM668A), centered on the SiO w = 1, J = 
1 — line (rest frequency 43122 MHz), and a broad band 
(50 MHz) centered on a line-free portion of the spectrum 
(offset in frequency from the maser by: ~350 MHz for 
AC952, -100 MHz for AC817 and AG622, ~40 MHz for 
AM668). Both frequency bands were observed in dual- 
circular polarizations. Absolute fiux density calibration 
was obtained from observations of 3C286 or 3C48 (de- 
pending on the epoch). Short (60 s) scans of the nearby 
(1.3°) QSO J0541-0541 (measured fiux density in the 
range 0.7-1.5 Jy, depending on the epoch) were alter- 
nated with 60 s scans of Source I to monitor amplitude 
gains, tropospheric and instrumental phase variations, 
and to determine the electronic phase offsets between the 
bands. ^ We then self-calibrated the narrow-band v = 1 
SiO maser data, applied the derived phase and ampli- 
tude corrections to the broadband continuum data, and 

^ NRAO is a facility of the National Science Foundation operated 
under cooperative agreement by Associated Universities, Inc. 

Project AM668A was conducted in A configuration plus the 
Pie Town VLBA antenna. 

3 In AM668A, the QSO 0501-019 was observed every 25 min- 
utes (relatively infrequently compared with the timescale of tropo- 
spheric and instrumental phase variations), so no absolute astrom- 
etry was available for that program. 
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TABLE 1 

Parameters of Observations. 



Program Date T Synthesized Beam (R = 0) RMS Synthesized Beam {R = 8) RMS 

(h) SmC) P.A(°) (mjy/bm) SmC) x e^C'); (mjy/bm) 



AMeeSA'^ 2000/11/10 6.6 0.041 x 0.028; -30 

AG622 d 2002/03/31 2.7 0.046 x 0.034; -4 

AC817 2006/04/15 0.9 0.051 x 0.030; 11 

AC952 2009/01/12 3.3 0.058 x 0.039; 3 



0.14 0.058 X 0.045; -20 0.13 

0.16 0.066 X 0.047; -17 0.14 

0.28 0.066 X 0.044; -5 0.24 

0.11 0.079 X 0.053; +4 0.11 



Note. — 

Approximate on-source integration time. 
^Synthesized beams correspond to images made using tlic AIPS task IMAGR witli a robust parameter R — and R — 8, 

respectively. 

'^Archival data from tire program AM668A liave been published bv lReid et all 120071) . 
'^Data from the program AG622 have been published bv lRodrfguez et al.l I I2005D . 



TABLE 2 

Parameters of the 7 mm sources in Orion BN/KL 



Sources a {J2000) S (J2000) et e„ Flux Density Angular Size 

™ (° ' ") (mas) (mjy) ^mC) x e™("); P.A.(°) 

BN 05 35 14.1094 -05 22 22.724 5 0.3 1 < 2 23 ± 2 0.095 ± 0.007 x 0.072 ± 0.008; 38 ± 10 

I 05 35 14.5141 -05 22 30.575 511<2 11 ±2 0.23 ± 0.01 x 0.12 ± 0.01; 142 ± 2 

n 05 35 14.3571 -05 22 32.719 5 4 <2 1.2 ±0.1 0.126 ± 0.018 x 0.064 ± 0.01; 6 ± 8 

H 05 35 14.5008 -05 22 38.691 5 2 <2 1.3 ±0.1 0.075 ± 0.015 x 0.053 ± 0.01; 4 ± 8 



Note. — The source names refer to Lonsdale et al. (1982), Garay et al. |198f|), and IChu rchwcll et al. (1987). The 
positions, flux densities, and sizes reported here are from the January 12 2009 epoch. Angular sizes were derived as 
described in the Text. The quoted uncertainties represent the dispersion among images with various weighting schemes. 
Individual contributions to the total error budget (~ 5 mas) on absolute positions are listed: et , en, , (see Appendix) . 



produced a high-quality map of the continuum emission. 
A detailed description o f this c alibration procedure can 
be found in iGoddi et all (|2009D . 

Using the Astronomical Image Processing System 
(AIPS) task IMAGR, we made images of the BN/KL re- 
gion with 4096x4096 pixels and cell size 0"005 to cover 
a 20" field. We produced maps with two different weight- 
ings of the (u,v) data, first setting the IMAGR "RO- 
BUST" weighting parameter to i? = 8 (equivalent to 
natural weighting) and then to R—0. We restored these 
image with a circular beams with FWHM 59 mas and 
47 mas, respectively. The sizes of these restoring beams 
are approximately the average of the geometric-mean 
sizes among different epochs (see Table [T]) . Since the 
synthesized beam size in the AC952 program was larger 
(^20%) than that in other programs, this choice resulted 
in a degradation of the resolution for other epochs. For 
this reason, we also produced images with round restor- 
ing beams of 40 mas (i? = weighting) and 50 mas 
[R — 8 weighting), equal to the geometric-mean size of 
the beams in projects AC817 and AG622, thus slightly 
"super-resolving" AC952. The "high-resolution" maps 
were used to analyze structural changes in the contin- 
uum emission from Source I (Sect. |4|), whilst the "low- 
resolution" images were used to fit positions and track 
proper motions of the continuum sources in BN/KL 
(Sect. [5|), since they minimize structural effects in lo- 
cating the peak of emission. 

3. RESULTS 

Four continuum sources were detected above a 6a lower 
limit of ^ 0.1 mJy bcam^^ in the most sensitive epoch 
(2009 January 12): BN. I. n. H (iLonsdale et all 119821: 
iGarav et al]|1987l: IChurchweh et al.lll987D . For each de- 



tected source, a two-dimensional ellipsoidal Gaussian 
model was fitted in a 100 x 100 pixel area around the peak 
emission, using the AIPS task JMFIT. For Source I, we 
used the AIPS task OMFIT to directly fit the visibilities 
using a disk model geometry, which is more appropri- 
ate for the elongated morphology of the emission than 
the Gaussian model fitted by JMFIT (see below). The 
positions, flux densities, and deconvolved angular sizes 
of the detected sources for the 2009 January 12 epoch 
are given m Table [2 The flux densities obtained for BN 
and Source I are in good agreeme nt with the values pre- 
viously reported in the literature (jMenten fc Reidl[l995l : 
iReid et all 120071 ) . source n was recent ly detected with 
the V LA for the first time at 7 mm bv iRodriguez et al] 
(|2009l ). who reports a flux density of 2 ± 0.2 mJy. This 
a factor ~2 higher than what we found in 2009 January 
12. The discrepancy might be due either to a different 
filtering of the total flux in the two observations (0"05 vs 
0"2 beam) or to intrinsic variability of source n, possi- 
bly due to non-therma l nature of the continuum emission 
(jMenten fc Reidlll995D . Based on a co mparison between 
7 mm and 3.6 cm observations, however. [Rodriguez et ahl 
(2009) estimated a spectral index of 0.2 for the contin- 
uum emission from source n, suggestive of marginally 
thick free-free thermal emission, as expected in an ion- 
ized outflow. Since the 3.6 cm emission is extended over 
0"6 (i.e., an order of magnitude more than the restoring 
beam in A-configuration), the variation in flux between 
measurements with different angular resolutions is most 
probably due to flltering of the total flux in the most 
extended configuration. 

Our images of the continuum emission from Source I 
at 43 GHz for the 2009 January 12 epoch are shown in 
Figure [U where different weightings and restoring beams 
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were used (59 mas, top panel; 47 mas, middle panel). 
The emission is composed of a compact component at 
the center and a component elongated NW-SE. In or- 
der to better demonstrate the structure of the elongated 
component, we have subtracted a 6 mJy (beam-sized) 
Gaussian component at the centroid of the emission. The 
resulting image (Fig. [U bottom panel) reveals a struc- 
ture with a radius of ~ 45 AU and a brightness of about 
4.5 mJy beam~^, elongated at a position angle (P.A.) of 
142 ±2°. It is barely resolved along the minor axis and 
has a deconvolved axial ratio of about 2:1. The radius 
of this elongated structure agrees to within uncertain- 
ties w it h previous measuremen ts reported bv lReid et al.l 
(|2007t ). iRodrfguez et al.l ()2009[ ) report a smaller size for 
Source I at 7 mm (0'.'09±0'.'01 x0'.'06±0'.'02). Their im- 
age, however, has lower sensitivity than the ones reported 
here, owing to shorter integration time and narrow band- 
width employed (line-free channels were used to map the 
continuum), hence the quoted size refers to the central 
compact component while the elongated component is 
not detected. 

Figure [2] shows a comparison among images of Source I 
at different epochs. Both size and position angle of the 
continuum emission appear to be remarkably constant 
over the last decade. 

Extreme care was used in computing the absolute as- 
trometry for the three new (non-archival) datasets re- 
ported here, the accuracy of which we quantify as ~5 mas 
(see Appendix). This accurate absolute astrometry en- 
abled us to track proper motions for BN and Source I over 
three epochs. Figure [3] shows the contour images of BN 
and Source I at different epochs. The displacement of the 
peaks in each image shows the absolute proper motions 
of the radio continuum emission for both sources. The 
consistency of morphology in the continuum emission at 
different epochs reinforces the robustness of the motion 
measurements. The proper motions have been calculated 
by performing an error-weighted linear least-squares fit 
of positions with time (Figure |4]) . The uncertainties of 
the absolute proper motions are the formal errors of the 
linear least-squares fits and depend mainly on the error 
of the absolute position of individual sources at differ- 
ent epochs. The derived absolute proper motions and 
associated uncertainties are reported in Table [H The 
measured proper motions are in the range 15-22 km s~^ , 
confirming the large proper mot ions of source s I and BN 
found previously at 8.4 GHz bv iGomez et al.l ()2008l ). In 
Table [3l we also report the absolute proper motions of 
both sources in the Orion Nebular Cluster (ONC) rest- 
frame and relative proper motion of BN with respect to 
Source I. The ONC rest-frame is obtained by subtracting 
the mean absolute proper motions of 35 radio sources lo- 
cated within a radius of 0.1 pc from the core of the ONC, 
/XqCOsJ = -|-0. 8ib0.2 mas yr~^; us = — 2.3±0.2 mas yr~^, 
PA=341±5° (|G6mez et al.l[2005h . 

Based on only two epochs, source n and H do not show 
significant proper motions in right ascension or declina- 
tion at a typical 3ct upper limit of 3 — 4 mas yr~^. 

4. MORPHOLOGICAL EVOLUTION OF THE RADIO 
CONTINUUM: AN IONIZED DISK AROUND SOURCE I 

The continuum emission from Source I shows an elon- 
gated structure with a position angle of 142° ± 3° and 
a size ^45 AU (Fig. [T]). The uncertainty quoted here 
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Fig. 1. — Continuum images of Source I in Orion BN/KL at 
43 GHz made with the VLA in the A configuration (2009 January 
12). The images have been produced with natural weighting and 
a round restoring beam of 59 mas {upper panel) and with a i? = 
weighting and a round restoring beam of 47 mas {middle and lower 
panel), respectively. In order to demonstrate the elongated struc- 
ture, we subtracted a (6 mJy, 50 mas size) Gaussian component 
at the centroid of the emission {lower panel). In all images the 
contour levels are at integer multiples of 0.4 mJy beam~^. The 
FWHM of the restoring beams are shown in the bottom right cor- 
ner of each panel. We interpret the northwest-southeast elongated 
emission, as coming from an ionized disk surrounding Source I (see 
also Reid et al. 2007). 

for the position angle represents the dispersion of values 
measured in several epochs, which is not significantly 
larger than the uncertainty in a single epoch (2°) re- 
ported in Table [2] A comparison among images at dif- 
ferent epochs does not show any evidence of changes in 
structure, size or position angle over 8 years (Fig. [2]). In 
addition, we found no evidence for significant variations 
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Fig. 2. — Continuum images of Source I in Orion BN/KL at 
43 GHz observed with the VLA in the A configuration at sev- 
eral epochs {upper panel). We subtracted the brightest component 
near the center in each epoch {lower panel). The images have been 
produced v^ith natural weighting and a circular restoring beam 
of 50 mas {bottom right corner). The contour levels are at inte- 
ger multiples of 0.4 mjy beam"! (2000), 0.5 mjy beam"! (2002), 
0.7 mJy beam-l (2006), 0.4 mjy beam"! (2009). The choice of 
contours is dictated by differing sensitivity and uv-coverage in dif- 
ferent experiments (Table [!}. There is no significant change in 
the morphology or position angle of the structure over ~8 years 
(142° ± 2°). 

in the total 7 mm continuum flux density of Source I 
over the last decade, with all values being consistent with 
11 ± 2 mJy. 

What is the nature of this ionized emission? The 
centimeter-to-millimeter wavelength spectrum of the en- 
tire source can be c haracterized as a pow er law with flux 
density S^, (x v ^■^ (|Beuther et al.ll2006[ ). The emission 
is optically thick up to 100 GHz, with a dust compo - 
nent dominating above 300 GHz ()Beuther et al.l [20061 ) . 
The radio emission has been inte rpreted eit her as an 
ioniz ed jet with a NW-SE axis (Ball v fc Zinne cker 2005; 
[Tan 2008a|) or a s an ionized disk with a NE-SW axis 
(|Reid et al]|2007( ). In the following, we will discuss three 
alternative interpretations of the radio continuum emis- 
sion from Source I: 1) an ionized jet, 2) an equatorial 
wind, and 3) an accreting disk. We will show that only 
the latter is consistent with our multi-epoch data. 

If the emission originates in an ionized jet, one would 
expect to detect structural variations owing to motions 



of individual jet components. In fact, when observed 
with high-angular resolution, radio jets do not show a 
smooth structure, but contain clumps. Using multi- 
epoch VLA observations, detection of proper motions 
as large as ^100-500 km s~^ have been reported in the 
components of radio continuu m jets associated with both 
low-mass YSOs (e.g., HHl-2 ^Rodriguez et al.ll200qD and 
high-mass YSOs (e.g., Cepheus-A - ICuriel et a'Ln 2006l: 
HH80-81 - iMarti et ed] 119951 ). Assuming 20 km s'^ as 
a lower limit for the expansion velocity of the jet asso- 
ciated with Source I (as estim ated from the SiO maser 
flow bv [Matthews et al.l[2010l ). displacements of --40 AU 
are expected on a temporal baseline of 9 years. Such 
displacements are comparable to the radius of the con- 
tinuum emission and are easily detectable with the linear 
resolution of our observations (^25 AU). In principle, no 
structural changes would be expected from an optically 
thin, homogeneous ionized region. However, even in the 
case of uniform density, thermal pressure would drive ex- 
pansion of the ionized region as a whole, a s seen for ion- 
ized winds and HII regions (e.g.. lAcord et al. 1998). As 
well, classic models of spherical winds ([Panagia fc Fellil 
Il975f l and ionized jets ()Revnoldslll986l ) predict the size 
of the ionized region to vary as a function of frequency 
as I/"" '', which would imply a factor of 3 change in size 
between 8.4 GHz and 43 GHz. Nonetheless, deconvolved 
sizes obtained for the radio continuum emission around 
Source I are similar at the two frequencies: 0''23 x 0"115 
at 43 GHz (this wo rk) and 0"19 x <0"15 at 8.4 GHz 
()G6mez et al.ll2008[ ). 

An alternative interpretation is that the radio emis- 
sion from Source I arises in a disk-like structure rather 
than in a jet, where the gas is not accreting but rather 
consists of an equatorial ionize d wind. Such an e quato- 
rial wind has been modeled by iDrew et al.l (|1998[ ) . The 
model predicts that the radiation pressure (mediated by 
line opacity) from a massive YSO accelerates material 
from the surface of the disk and blows it away mainly in 
the equatorial plane. The velocities obtained for the gas 
are of the order of a few hundreds km . Also in this 
scenario, large changes in morphology are expected over 
short time intervals (i.e., a few years) owing to fast mo- 
tions of clumps in the wind, a s observed in the luminous 
YSO S140 IRSl (|Hoarel 120061 ). In addition to fast mo- 
tions, rapid changes in the wind structure are expected 
as well, owing to small changes in illumination by the 
ionizing YSO and to the inherently unstable nature of 
the radiation-driven winds (Owocki et al. 1988). Hence, 
we conclude that our multi-epoch data are not consistent 
with the picture of a radiation-driven disk-wind model. 

Finally, we consider the case that the elo ngated struc- 
ture traces ionized accreting gas in a disk. iKetol (|2003l ) 
has shown that, even when the inner portion of a disk 
is fully ionized, accretion of (ionized and neutral) mate- 
rial can continue onto a massive YSO. In this spherically 
symmetric model, the gravitational attraction of the star 
exceeds the outward force of the (radiation and thermal) 
pressure inside a critical radius = GM/2cs^ where G 
is the gravitational constant, M is the mass of the YSO, 
and Cs is the sound speed in the ionized material. In the 
hypothesis of spherical acc retion, by s caling the stellar 
parameters in Table 1 from lKetol (2003) to a 10 Mq star, 
ionized accretion can proceed inside of the critical ra- 
dius '^25 AU with an accretion rate ^ 10^^ yr^^. 
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This critical radius is reasonably consistent with the disk 
size observed at 43 GHz in Source I. In Sect. [6l we re- 
port evidence for Source I being a binary of total mass 
of the order of 20 Mq . In the most probable case of an 
equal-mass binary, this would slightly increase the crit- 
ical radius above, estimated in the assumption of one 
10 M0 single star. We note that similar values for the 
critical radius and accretion rates are obtained for models 
of ionized ac cretion in th e presence of significant angular 
momentum (|Ketol 12007*) . as required by the disk geom- 
etry in Source I. We also note that an accretion rate of 
10~^ Mq yr~^and a central stellar mass of 10 to 20 Mq is 
consistent with a luminosity of a few times 10"* Lq for 
Source I. 

The interpretation of the ionized elongated structure 
as a disk versus a jet is relevant to a key longstanding 
question regarding Orion BN/KL: what are the drivers 
of the region's outflows? In our interpretation, the elon- 
gated radio continuum structure traces a nearly edge- 
on disk, with a spin axis aligned NE-SW, parallel to 
the low- velocity outflow (see § [T|). This is consistent 
with the model for Source I based on high - resolution 
monitoring of several SiO maser transitions pCim et al.l 
l2008HMatthews et al.|[2010l: iGreenhill et al.ll2004al and in 
prep). Greenhill et al. (in prep.) estimated a mass-loss 
rate of the NE-SW outflow to be 3 x 10"^ Mq yr-\ 
assuming a gas density of 10^ cm~^ (required for the 
excitation of w = SiO masers), a typical velocity of 
18 km s~^ for the outflow, and a radius of 100 AU. Based 
on the models of iKetd (|2007( ). we estimate an accretion 
rate ~ 10~^ M© yr~^. A ratio of 0(0. 1) for the mass- 
loss rate to the accretion rate is typical of values inferred 
from observations of l ow-m ass pre-main-sequence stars 
(e.g., iBontemps et al.l 19961) and of values predicted in 



MHD ejection models fehu et al.ll2000HKonigl fc Pudritzj 
|2000( ). This might indicate that Source I is still in an 

active accretion stage. We caution that the mass-loss 
rate of the outflow from Source I estimated from the SiO 
masers is accurate only to within an order of magnitude 
(Greenhill et al., in prep.). Nevertheless, the preponder- 
ance of evidence strongly supports a picture where the 
high-mass YSO Source I is powering a disk/outflow sys- 
tem along a NE-SW axis. 



5. PROPER MOTIONS: A CLOSE PASSAGE BETWEEN 
SOURCE I AND THE BN OBJECT 

We observe BN and Source I moving with high 
speeds in opposite direction with respect to one an- 
other. This result confirms previous i ndependent proper 
motion measurements at 8.4 GHz (iGomez et al.l [20081 ) 
and at multiple freque ncies ranging from 8.4 to 43 GHz 
(jRodriguez et al.|[2005() . Figure [5] shows maps of the ab- 
solute proper motions of both sources in the ONC rest- 
frame (upper panel) and relative proper motion of BN 
with respect to Source I (lower panel). 

The more accurate relative proper motions show that 
Source I falls within the error cone of the BN past motion 
(Fig. [SI lower panel). Based on our accurate relative 
astrometry, we can estimate the minimum separation BN 
and Source I had in the past in the plane of the sky 
by extrapolating the velocity vector back ward in time 
and assuming no acceleration. Following iGomez et al.l 




KIGHT ASCKNSION (J^OGO) 




RIGHT ASCENSION (J3000) 



Fig. 3.— Contours for the 2002 (blue), 2006 (red), and 2009 
(black) epochs of radio source I (upper panel) and the BN object 
(lower panel). Contours are (1, 2, 3, 5, 7, 9, 11) X 0.5 (2002), 0.7 
(2006) and 0.4 (2009) mjy beam"! (upper panel), and (1, 3, 6, 9, 
12, 15, 18, 20) X 0.7 (2002), 0.9 (2006) and 0.7 (2009) mJy beam-i 
(lower panel). The images have been restored with a circular beam 
with FWHM 59 mas (bottom right corner). The displacement of 
the peaks in each image shows the absolute proper motions of the 
radio c ontinuum emission of Source I and BN. 

(|2008f ). the minimum separation is 



F2OO2.25M1/ ^ 2/2002. 25Ma: I 



which occurs at an epoch t„ 
ri„in(years) = 2002.25 - 



given by 



|x2002.25Aix + y2002.25My I 



where both quantities are expressed as a function of 
the relative position for program AG622 (0:2002 25 = 
-5.974" ± 0.002", 2/2002.25 = 7.752" ± 0.002") and rel- 
ative proper motions of BN with respect to Source I 
= -0.0107" ± 0.0003"yr-i, and ^ly = 0.0142" ± 
0.0003"yr~-'^). By substituting these values, we obtain: 
S^i^ =0"109±0"175 and T^in = 1453 ± 9 years. Hence, 
our new measurements indicate that, about 560 years 
ago, BN and Source I were as close as 50 ± 100 AU in 
projection on the plane of the sky. 
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TABLE 3 

Proper, Motions of Radio Source I and the BN object 



Sources tJ-s Mtot Vy V^od P-A. 

(mas yr"-"-) (mas yr"-"-) (mas yr~-^) (km s"-"- ) (km s^-^ ) (km s"-"- ) (deg) 

Absolute proper motions 

I 6.3 ±1.1 -4.2 ±1.1 7.6 ±1.1 12.3 ±2.1 -8.3 ±2.1 14.9 ±2.1 123.9 ± 8.2 

BN -3.5 ±1.1 10.4 ±1.1 11.0 ±1.1 -6.9 ±2.1 20.4 ±2.1 21.5 ± 2.1 341.3 ± 5.5 

Proper motions relative to I 

BN -10.7 ±0.3 14.2 ±0.3 17.8 ± 0.3 -21.0 ±0.6 27.9 ± 0.6 34.9 ± 0.6 323.0 ± 1.0 

Proper motions in the ONC rest-frame 

I 5.5 ±1.1 -1.9 ±1.1 5.8 ±1.1 10.8 ±2.1 -3.8 ± 2.2 11.5 ±2.1 109.4 ± 10.9 

BN -4.3 ±1.1 12.7 ±1.1 13.4 ±1.1 -8.5 ±2.1 25.0 ±2.1 26.4 ±2.1 341.2 ±4.6 



2004 2006 
Epoch (years) 



2004 2006 200B 
Epoch (years) 



Fig. 4. — Measured absolute proper motions in RA and DEC of 
the radio continuum centroid of Source I and BN, relative to po- 
sition a(J2000) = 05'"35™14=.5116, 5(J2000) = -05°22'30!'536. 
In each panel, the dotted line shows the proper motion calculated 
by the (variance-weighted) linear least-squares fit of the positional 
offsets with time. 



Although we do not have position information along 
the line-of-sight that confirms the close passage in 3- 
dimensions (proper motion analysis gives only a 2- 
dimension description of the dynamics), the line of 
sight velocities of BN and Source I are also con- 
sistent with a dynamical encounter. Based on the 
measured radial velocitie s of the ambient mol ecular 
cloud (yLSR=9 km s-i : IGenzel fc StutikH HOSl) . BN 
(^LSR=20 km , or 11 km s~^with respect to the 
ambient cloud: [Rodriguez et al.l I2OO90 . and Source I 
(^LSR=5 km s~^ , or —4 km s~^ wit h respect to the am- 
bient cloud: [Matthews et al.ll2010[ ). the ratio of radial 
velocities between BN and Source I is ^^2.7, which is con- 
sistent with the value of 2.3 determined from the (ONC- 
frame) velocities in the plane of the sky. Hence, 5 to 6 
observables are consistent with the hypothesis of a close 
passage between Source I and BN. We also note that 
the closest approach separation of 50 AU is comparable 
to the radius of the ionized disk observed today around 
Source I (Sect. 0]). As this is compact for a high-mass 
YSO accretion disk, we suggest that it may have been 
truncated in the close encounter with BN. We note also 
that the ionized envelope surrounding BN has a compa- 
rable radius of ^ 40 AU (Table [2]), consistent with the 



tr uncation hyp o thesis . 

I Gomez et al.l ()2008[ ) measured a large proper motion 
for source n as well (/i = 13 ± 1 mas yr~^, v — 
25 ± 2 km s^i , P.A.= 180° ± 4°) using position mea- 
surements at 8.4 GHz and proposed a close passage with 
BN and Source I. This result however disagrees with our 
measurements. Using two epochs in which we achieved 
5 — 6(7 detections at 7 mm wavelength, we measured 
a much smaller proper motion of 3.7 ± 1.6 mas yr^^ 
and PA= 26 ± 22 in the ONC rest-frame. This in- 
dicates that source n is moving either very slowly or 
not at all (our measurement is consistent with within 
~ 2a = 3.2 mas yr~^). In addition, the proper motion 
vector points towards and not away from the interaction 
center, as it would be required if source n were truly 
ejected by a close passage with Source I and BN. As a 
final point , the proper motion of source n at 8.4 GHz re- 
ported bv lGomez et al.l ()2008[ ) should be viewed with a 
degree of caution. The emission at 8.4 GHz from source n 
is resolved into a bipolar or elongated structure, depend- 
ing on epoch (size~0"6), and the systematic error in 
the motion may arise from internal variability. Track- 
ing of northe r n and southern components of source n by 
iGomez et al.) (I2n08h may be possible, but a longer time 
baseline is requ ired for certainty. Of some concern, the 
dataset from iGomez et al.l ()2008f ) is not homogeneous 
in terms of observing set-up among epochs, and hence 
changes in uv-coverage could mimic cha nge in struc- 
ture/p osition. Finally, in their Figure 2, IGomez et alj 
()2008[ ) show that the displacement of source n in the 
plane of the sky from the cluster origin is approximately 
half that of BN. However, this is inconsistent with the 
two objects having comparable proper motions and start- 
ing their outward motion from the same point of origin 
and epoch (assuming linearity and constancy of motion). 
Further observation of source n is required to confirm the 
proper motion measured at 8.4 GHz and hence the dy- 
namical interaction with BN and Source I. 

In the following section, we will reconstruct the dy- 
namical history of BN/KL based on the proper motions 
measured for Source I and BN. Section [6.3.2l will discuss 
the possible dynamical role of source n, in the assump- 
tion that the 8.4 GHz proper motion measurement of 
Source n derived by IGomez eTaJ] (f2008h is correct. 

6. DYNAMICAL HISTORY OF BN/KL 

Observations have increasingly shown that stars are 
almost never born in isolation but favor formation in 
clusters, where interactions can occur among members 
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Fig. 5. — Proper motions of BN and Source I measured from our 
multi-epoch datasets in the ONC rest-frame {upper panel). Proper 
motion of BN relative to Source I {lower panel). The arrows in- 
dicate proper motion direction and displacement for 200 yr. The 
dashed lines indicate errors in proper motion position angle, and 
they track backwards showing past positions of source. About 
560 years ago, BN, Source I, and several IR sources within the 
BN/KL region were located within a few arcseconds of each other. 
Our measurements are consistent with a null proper motion for 
source n. 



of a collapsing protostellar group. Interactions with 
sibling stars play an even m ore important role in th e 
formation of high- mass stars ()Ballv fc Zinneckeill2005fl . 
YSOs are also surrounded by circumstellar disks and en- 
velopes, which can play a role in this context by increas- 
ing the interaction radius amo ng members of the co l- 
lapsing protostellar cluster (e.g. iMoeckel fc Ballvll2006l) . 
The anamolously large relative motions between BN and 
Source I suggest that interactions have also played an 
important role in Orion BN/KL. In the following sec- 
tions, we investigate four scenarios to account for the 
observed dynamics of the region: 1) the motions are in- 
dependent and BN is a runaway star from the Trapez- 
ium cluster fSect l6.1) ): 2) a close passage between the two 
protostars resulted in a star-disk interaction (Sect 16. 2p ; 
3) the motions are the result of a multiple (> 3) body 
interaction within a cluster of individual YSOs, where 
a tight massive binary is formed (Source I) that carries 
away the binding energy of the system and BN is ejected 
(Sect l6.3|) : 4) a close passage between BN and Source I, a 
pre-existing massive binary, resulted in the acceleration 
of both stars and the hardening of the original binary 
(Sect EH). 



6.1. Scenario 1. 



Dynamical ejection of BN from 6^ 
C 



Ori 



We consider first the possibility t hat the motions o f 
BN and Source I are independent. [S3 (|2004 l2008bf ) 
proposed that BN was ejected from the Trapezium by 
dynamical interaction involving the 9^ Ori C binary sys- 
tem. 

In a triple system usually the lightest member is 
ejected, sometimes as a runaway star, while an eccen- 
tric binary is formed between the remaining stars. The 
formation of a tight binary is required if the system was 
originally bound, as a consequence of energy conserva- 
tion. Usually, both components (star and close binary) 
leave their parental envelope, with velocities inversely 
proportional to their mas ses, as a consequence of mo- 
mentum conservation Ce.g. lReipurth| [2000). 

6*^ Ori C seems to have the physical properties pre- 
dicted for such a dynamical event: an eccentric, mas- 
sive binary with p rimary and secon dary masses greater 
than that of BN (Kraus et al.ll2007l) : a total orbital en- 
ergy {Etot ~ 3 X 10"*^ erg) greater than the kinetic en- 
ergy of BN (10'*^ erg); a proper motion roughly in the 
opposite direction to the motion of BN and with ap- 
proximately the predicted magnitude (for the dynami- 
cal mass of BN). However, our absolute proper motion 
measurement is inconsistent with a close encounter be- 
tween BN and 9^ Ori C in the past. Figure |6] shows the 
radio-ONC-frame^ proper motion of BN, = —4.3 ± 
1.1 mas yr-\ ^Xy = -f 12.7 ± 1.1 mas yr~i (Table E]), 
and the optical-ONC-frame^ proper motion of 9^ Ori C, 
fix = +1.4 ± 0.17 mas yr iiy = —1.8 ± 0.16 mas yr ~^ 



* The radio-ONC-frame is defined by the mean absolute proper 
motions of 35 radio sources located within a radius of 0.1 pc from 
the core of the ONC, fiaCos5 = -1-0.8 ifc 0.2 mas y r~^ ; fig = —2.3 it 
0.2 mas yr"!, PA=341 ± 5° lIGomez et al.|[2005r ). 

^ The optical-ONC- frame is defined by the mean absolute proper 
motions of 73 stars brighter than V ~ 12.5 located within a radius 
of one-half degree of the Trapezium, correspond ing to an average 
prope r motion dispersion of 0.70it0.06 mas yr~^ dvan Altena et all 
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(|van Altena et allll988[ ). The plot shows that past tra- 
jectories of BN and 9^ Ori C in fact intersect far from 
the Trapezium cluster. By repeating the same analy- 
sis we made for the relative motions of BN with re- 
spect to Source I (this time referring positions and ve- 
locities to 6*1 Ori C^), we obtain xbn = -35.4" ± 0.3", 
ysN = 60.0" ±0.3" , fi^ = -0.0057" ±0.0011" yr-\ and 
Hy = 0.0145" ±0.0011" yr-i (P.A.=339 ± 4°) for the rel- 
ative proper motions between BN and 9^ Ori C (today). 
This yields S^in = 10.9" ± 4.8" and T^in = -2403 ± 312 
years. Hence, our measurements indicate that probably 
BN and 9'^ Ori C did not get closer than -5000 AU in 
the past. This large separation is much greater than the 
gravitational radius that corresponds to the inferred es- 
cape speed of BN, which argues against a close encounter. 
Nevertheless, the accuracy of our measurement (at a 2.3cr 
level) does not conclusively exclude it. We note that the 
un certainties her e are a factor of 5 larger than quoted 
by iGomez et al.l ()2008[ ). whose analysis used the rela- 
tive motion of BN with respect to Source I to achieve a 
greater accuracy. However, since Source I has a large ab- 
solute proper motion (see Table |3]), using the BN motion 
relative to Source I is undesirable. 

Apart from the measurements reported in this paper, 
there is a number of pieces of evidence that are inconsis- 
tent with an ejection of BN from 9^ Ori C. First, 9^ Ori C 
and BN have similar absolute radial velo cities mea- 
sured in the ONC rest- frame, —15 km s ~^ (|Stahl et al.l 
|2008[) and 11 km s~^ (|Rodriguez et al.l I2009D . respec- 
tively, which is inconsistent with momentum conserva- 
tion, since BN is approximately five times less massive 
than 9^ Ori C. Second, a close passage of BN and 9^ Ori C 
would explain the large proper motion of only BN but 
not that of Source I. We note that in order to reach a 
velocity dispersion of 15 km s~^ , a stellar cluster on the 
order of 1000 AU in extent would have to reach 250 M0, 
which is greater than the mass concentrations believed 
to exist in the part of BN/KL local to Source I. Third, 
this scenario would not explain the close passage (~0"1) 
of BN with Source I, since it would be a notable coin- 
cidence for a high-velocity YSO ejected from a massive 
binary to pass so close to another high-mass YSO. Fi- 
nally, BN exhibits characteristics of a star much younger 
than any of t he Tra pezium members, which according to 
iHiUenbrandl (|l99l are a few million years old. While 
none of the Trapezium stars are surrounded by obvious 
circu mstellar material, BN is surrounded by a molecular 
disk (|Scoville et al.lll983t iBeuther et al.ll2010D and by a 
hypercompact HIT region with a ~20 AU radius, imply- 
ing that this star is dragging along a significant amount 
(> 10 Mp)) of dense circum stellar gas (Scoville et al.l 
119831: iBallv fc Zinneckedl200l . 

6.2. Scenario 2: Two-body interaction: star-disk 
encounter with BN 

We now consider the case where BN and Source I were 
initially unbound and approached to within ^^50 AU 
(Sect. [SJ. Source I is surrounded by an ionized disk, 
indicating a massive (possibly still accreting) protostar 
(Sect. 14]) . A close passage by a massive star (BN) should 

^ For 9^ Ori C we took the absolute position from the COUP X- 
ray sur vey: 05''35"'16.478y , 05°23'22.844", known to an accuracy 
of 0!'3 jGetman et al.|[2005l) . 




5''35'°17° 16' 15' 14' 

RIGHT ASCENSION (jaOOO) 



Fig. 6. — Absolute proper motions of BN (radio, this work) 
and 9^ Ori C (optical, Ivan Altena et all [198 81 measured in the 
ONC-rest-frame. The arrows indicate proper motion direction and 
displacement on the plane of the sky over 1000 years. The dash- 
dot line indicates motion backward over 4400 years, assuming it 
is linear. The dashed lines indicate 1 a uncertainties in position 
angles. The minimum separation on the sky between BN and 9^ 
Ori C was 10.9" ± 4.8". 

produce tides in the disk a nd enhan c e acc retion (e.g. 
lOstrikeHfl99l . In particular, |Pfalzner| ([20061) has shown 
that disks around massive YSOs located in clusters can 
experience angular momentum losses of up to — 50%-95% 
(with consequent increases in accretion), as a result of 
gravitational interaction during stellar encounters; this 
is termed cluster- assisted accretion. 

Although star-disk interaction at the center of clus- 
ters can affect accretion, it cannot itself lead to stellar 
ejection. In fact, a two-body interaction does not allow 
efficient transfer of momentum and/or energy between 
the two objects. On the contrary, simulations show that 
disk-assisted capture i ncreases binary formation rates 
(jMoeckel fc Ballvl2007D . In order to account for the large 
(positive) total energy of the BN-I system, the formation 
of a tight binary that develops a large negative-binding 
energy is required, a consequence of momentum and en- 
ergy conservation. 

In conclusion, in order to end up with stellar ejections, 
the original system must have contained three or more 
objects, and one of the successor systems must now be a 
tight binary. 
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TABLE 4 

Summary of dynamical scenarios for BN/KL and N-body simulation outcomes. 



Scenario 






Ejections'^ 


Smin Passage"* Sect. 




(Mq) 


(Mo) 


(%) 


(AU) 



3-body: * * * =^ * (* *) 
BN, I (new binary) 

4-body: 

BN, I, n (new binary)° 



10,10,10 



10,10,5,5 
10,10,9,1 



10+10 



5+5 
9+1 



11-body: llx(*) ^ 9x(*) {* *) 0.5(x2),l(x2),2(x2) 10+10 
IRcx,BN,n,I (new binary) 3(x l),5(x 1),10( x3) 10+5 



3-body: * (* *) ^ * (* *) 
BN(single), I (old binary) 



10,10+10 



10+10 



0.1 



< 0.1 

< 0.1 

0.2 
0.2 

38 



.3.2 



.3.2 



5 
4 

16 



6.3.3 



[631 



Note. — 

''Stellar masses in the N-body simulations. 

''Masses of binary components that remain after ejection of stars from the system. 
"^Fraction of cases leading to the ejection of a massive (~10 Mq) object. 

Average mimimum passage (and binary separation) for cases that include an ejection. 
°It assumes that n is a binary with a 5+5 or 9+1 total mass, source I and BN two single 10 stars. 



6.3. Scenario 3: Three and more-body interactions: 
single stars 

In this Section, we describe different scenarios that can 
result in the high proper motion of multiple stars: a triple 
system (Sect. 16.3.ip . a quadruple system (Sect. 16.3.2p . 
and a stellar cluster (Sect. 16.3. 3p . In order to test 
different scenarios and quantify statistical probabilities, 
we r an N-body sim ulations using the NBODYO code 
from lAarsethl ()1985f ).^ We numerically integrated sev- 
eral thousand cases of cluster decays with different initial 
numbers of bodies. Table |4] summarizes the key aspects 
and likelihood of different scenarios, which will be de- 
scribed in detail in the following subsections. We note 
that each scenario involves the ejection of individual stars 
and concomitant formation of a tight binary. The pres- 
ence of a disk around Source I poses a challenge for some 
arrangements (Sect. 16.41) . The existence of a binary be- 
fore the encounter provides an evolutionary path that 
may agree best with the available data (see Sect. I6.5|) . 

6.3.1. Decay of a triple system and formation of a new 

binary 

In the simplest dynamical scenario of a triple system, 
we assume that Source I, which is apparently moving 
slower than BN, forms an eccentric massive binary, and 
lea ves the system al o ng wi th B N, as previously prop osed 
by iRodriguez et al.l (|2005[ ) and iGomez et al.l ()2008[ ). In 
order to explain the observed velocities, we need to es- 
timate the mass of both objects. Based on a luminosity 
of 25 00 Lq for BN from mid-IR imaging (|Gezari et alJ 
|1998( ) and an ionizing photon rate of ~ 4 x 10^^ from 
optic al ly thin emission a t v = 216 GHz (jPlambeck et al.l 
119951 ). IRodriguez et al.l (|2005f ) report a mass Mbn = 
10 ± 2 M0. In the following, we will assume a fiducial 
value of 10 M0. Constraining the mass of Source I based 
on an infrared or bolometric luminosity is hindered by the 
very high extinction toward the source (jGreenhill et al.l 

NBODYO is a direct N-body integrator which for each particle 
computes the force due to all other N-1 particles. Each particle is 
followed with its own integration step, an essential feature when 
the dynamical timescales of different particles in the simulation 
could vary significantly. 



l2004bD . iMatthews et al.l (|2010[ ) estimate a dynamical 
mass of >7 M0 from SiO masers, whereas 'Rei d et all 
j2007l ) report a mass ~ IOMq based on modeling of 
proton-electron Bremsstrahlung at 7 mm. However, 
from conservation of linear momentum for the BN-I sys- 
tem, one would estimate a mass: Mj — 2.3 x Mbn ^ 
23 M0 for Mbn = 10 Mq. This result requires explana- 
tion of systematic errors i n estimation of the dyna mical 
mass from SiO masers bv IMatthews et al.l (|2010l ). We 
discuss implications in Sect [71 

Given mass estimates for the stars, the orbital sepa- 
ration of the binary can be estimated from conservation 
of mechanical energy. By assuming that the total ki- 
netic energy of the system today, ^{Mjvj + MBNv'^j^) ^ 
lO''^ erg, is approximately equal to the final binding en- 
ergy of the binary, {GMjpMjs)/2a, where Mjp and M/^ 
are the masses of the primary and secondary stars, re- 
spectively, and a is the binary semi-major axis, we ob- 
tain (a/AU) = 40mf(l - mf)(MBN/10 M©), where m/ 
is the primary's mass fraction in the binary (to/ = 
Mip/{Mip + Mis)). In the assumption of an equal-mass 
binary (to/ = 0.5), the semi-major axis is 10 AU for 
Mbn = 10 M0, corresponding to a period of '^7 yr. 
We note that an unequal-mass binary with the same to- 
tal mass would imply an even smaller orbital separation, 
but it would be dynamically less favorable for ejection 
(see discussion in following sections). 

Hence, a ^ 20 M© equal-mass binary with a semi- 
major orbital axis < 10 AU will have a binding energy 
sufficient to compensate the positive kinetic energy of the 
BN-I-Source I system. Direct observation of such a tight 
binary may be beyond what can be achieved with current 
data. Size scales of 10 AU cannot be resolved with the 
VLA, and no maser emission has been detected with the 
VLBA at s uch small radii (resulti ng in the "dark band" 
reported bv IMatthews et ani201Clf ). However, highly ec- 
centric orbits may be anticipated to generate local fluctu- 
ations in disk heating over just one quarter orbital period 
(~ 2 years), and these might be recognizable over time. 

The simple analysis above allows estimation of phys- 
ical parameters of the binary using basic physics (mass 
from linear momentum conservation and orbital sepa- 
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ration from mechanical energy conservation). In order 
to quantify the probability of such a triple-body inter- 
action, we simulated a system composed of three single 
stars of 10 Mq mass randomly distributed inside a sphere 
of 1000 AU diameter. This corresponds approximately 
to the area subtended by the error cones of proper mo- 
tions of BN and I (Fig. 5, top panel). We ran 1000 
simulations for durations of 800 years, adopting a veloc- 
ity dispersion of 0.4 km s^^ , corresponding to a thermal 
velocity in a molecular cloud at T = 10 K. We found that 
only 0.1% of configurations (1 out of 1000) gave ejections 
with speed higher than 25 km s~^ . In the case with ejec- 
tion, a 10-1-10 M0 binary forms with orbital separation 
~ 7 AU (Table m). The single 10 Mq star is ejected with 
a velocity 30 km s~^ , while the massive binary leaves 
the interaction area with a velocity ^ 15 km , as ex- 
pected from conservation of linear momentum. These 
simple simulations show that the BN-I system is un- 
likely (probability 0.1%) to arise from chance encounters 
among three single stars. 

6.3.2. Decay of a quadruple system and formation of a new 

binary 

Source I and the BN object may not be the only mas- 
sive objects in BN/KL. Radi o source n, a YSO associ- 
ated with a bi polar iet fiMenten fc ReidI I1995D and an 
accretion disk (jGreenhill et al.l l2004b[ ). has been pro- 
posed as a possibly significant source of luminosity in 
BN/K L (iMenten fc ReidI 119951) . Using MIR measure- 
ments, iGreenhill et al. (|2004br i estimated a luminosity 
of 2 X 10'^ Lq, corresponding to a stellar mass ~ 8 Mq. 

Based on posit ion measurements at 8.4 GHz, 
iGomez et al.l (|2008D proposed for source n an origin 
in the cluster decay along with BN and Source I. In 
the following, we investigate this possibility by assum- 
ing a proper motion vxn = — 1.6 km s~^ a n d vun = 
-21.4 km s~i , as derived by iGomez et all ()2008f ) at 
8.4 GHz in the ONC rest-frame (but see discussion in 
Section [5] for arguments as to w hy this measurement is 
questionable). For M„ = 8 Mq (jGreen hill et al.ir2004b[) 
and Mbn = 10 Mq (Sect. I6.3.1|) . conservation of linear 
momentum (along right ascension and declination) gives 
a mass Mi ~ 8 Mq, much lower than that estimated for 
the triple-system decay (Sect. [673TT|) . In this scenario, at 
least three YSOs (BN, I, and n) with comparable masses 
(~8-10 Mq) would be running away from each other with 
high velocities, after a dynamical interaction occurred 
about 500 years ago. 

However, formation of a binary is required to enable 
ejection of stars from a cluster. It is unlikely to be 
Source I because an equal-mass binary of 10 Mq to- 
tal would not produce ionizing radiation, while a binary 
with very different masses for primary and secondary is 
an unlikely outcome for an N-body interaction (see be- 
low). Hence, in the following analysis, we assume that 
n is the binary. A couple of arguments support this hy- 
pothesis. First, source n has been identified as a hard 
X-ray source, and since a B-type star is unlikely to be 
active, n might be an unequal-mass binary with one low- 
mass member. Second, based on J and H magnitudes, 
iLagrange et al.l ()2004[ ) stated that the SED of source n 
cannot be modeled as a high-mass star with a single ex- 
tinction value. An additional source of IR luminosity 
must be added (e.g., circumstellar dust or stellar com- 



panion). Although J and H band observations have not 
revealed bright companion IR sources with separations 
down to ~20 AU, a tighter binary remains a possibility. 

In order to estimate the probability that three ~ 
10 Mq stars are ejected after interaction, we simulated 
two 4-body interactions that could generate a 10 Mq bi- 
nary: 1) two 10 Mq and two 5 Mq stars and 2) two 
10 Mq, one 9 Mq, and one 1 Mq stars. Initial con- 
ditions were set as described in Sect. 16.3.11 Among 
1000 simulations, not a single event ends with ejection 
of a 10 Mq star and formation of a binary 5-1-5 Mq or 
9-1-1 Mq. In system 1) we observed ejections in 0.8% of 
cases with masses 5 Mq (0.7% of cases) and 10 Mq (0.1% 
of cases). The new binaries formed in those cases had 
total masses of 10-1-5 Mq and lO-flO Mq in four cases 
(0.4% of total) each. Orbital separations varied in the 
range 2-11 AU (average 9 AU) and 11-23 AU (average 
14 AU), for the 10-1-5 Mq and 10-1-10 Mq binaries, re- 
spectively. In system 2) we observed ejections in 2.8% of 
cases but only one with mass 9 Mq (0.1% of cases), when 
a binary 10-1-10 Mq formed with an orbital separation 
of 10 AU (the low mass member was ejected in almost 
all cases). Even assuming source n was a binary before 
the encounter (with component masses of 5-1-5 Mq or 
9-1-1 Mq), an exchange interaction would eject the low- 
mass companion and form a quasi-equal mass binary of 
^ 15 — 20 Mq. Hence, unless source n is about twice 
as luminous as believed, the star is unlikely (probabil- 
ity < 0.1%) to be the binary that enables ejections from 
4-body interaction (Table I5- 

The conclusion is co nsistent with other pieces of evi- 
dence. IGreenhill et al.l ()2004bD detected a disk around n 
with a radius ~ 170 AU, larger than the separation esti- 
mated for the BN-I encounter (50 ± 100 AU), indicating 
that either source n did not take part in the encounter 
at all or that the passage occurred at large separation 
(>200 AU). In this last case, source n could not have ex- 
changed momentum/energy e fficiently, making this hy- 
pothesis unlikely. In addition, ILagrange et al.l (|2004f ) re- 
port detection of a brown dwarf separated by ^--^300 AU 
from source n. If this is indeed a bound system, then it is 
unlikely to have survived an earlier dynamical interaction 
with BN and Source I. 

We conclude that, even assuming that the 8.4 GHz 
proper motion me asurement of source n derived by 
IGomez et all (|2008f) is correct, source n is unlikely to 
have played a major role in the dynamical interaction 
between Source I and BN. 

6.3.3. Decay of a large-N cluster 

To exp lain the observed motio n of BN with respect to 
Source I. lRodriguez et al.] (j2005l ) proposed an interaction 
among several members of a collapsing protostellar clus- 
ter. In addition to the radio sources studied here, many 
individual peaks of thermal IR (7 — 24/Ltm) emission have 
been identified in the BN/ KL region from high angular 
resolution IR studies (IGezari et al.lll998l : IGreenhill et al.l 
l2004bt IShuping et al.l 120041 ). So me of these may be em- 
bedded, self-luminou s protostars (|Greenhill et al.ll2004bl : 
IShuping et al.l |2004|) , suggesting a protostellar density 
even higher than that of the Trapezium cluster. For ex- 
ample, the 5 knots in IRc2 plus Source I would constitute 
a cluster core of 10*" — 10* pc^"^ over ^ 10"^ AU wi thin the 
larger (- 20") IR cluster (jGreenhill et al.ll2004bD . Apart 
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from IRc2, also IRc7 and IRc4 (and possibly IRc3) ap- 
pear to be clustered around the putative expansion center 
within a few arcseconds (Fig. [SJ. 

We propose that the original bound cluster may dis- 
integrate resulting in the expansion of multiple objects 
with different velocities and the recoil of one (or more) 
of the protostars as tight binaries. This scenario would 
explain why the putative center of expansion is presently 
devoid of infrared sources (see Fig. [S|). We explicitly 
note that, owing to the lack of astrometric IR observa- 
tions, the possibility that the members of the IR cluster 
in BN/KL are participating in the expansion measured 
for Source I and BN, remains speculative. 

In order to test the present scenario, we ran simulations 
of a protostellar cluster with eleven members, adopting 
masses of 0.5 (2 stars), 1 (2 stars), 2 (2 stars), 3 (1 star), 
5 (1 star), and 10 (3 stars) Mq (total mass of 45 Mq), 
randomly distributed inside a sphere of 1000 AU diame- 
ter. The implied stellar density is 1.6 x 10^ stars pc~^, 
comparable to that found in the Orion Tra pezium and 
the Cepheus A IIW2 stellar com plex (e.g. iComito et al.l 
l2007t iJimenez-Serra et aLll2009[ ). The choice of masses 
reflects plausible masses for the objects present in the 
region, based on radio and IR data. The three stars 
with 10 M0 are meant represent BN, I, and n, while the 
low-to-intermediate mass objects represent the IR knots 
associated with IRc2, IRc3, IRc4, and IRc7. 

We ran 1000 simulations with durations of 800 years, 
adopting a velocity dispersion of 0.4 km s^^ . In a sig- 
nificant fraction of configurations (17%), stars escaped 
from the cluster with velocities greater than that of BN 
(25 km s^^); we define these as runaway stars. The 
runaway stars had a maximum mass of 10 Mq in 0.4% 
of cases. High-mass runaways always required the for- 
mation of a binary with 10-1-5 Mq (0.2% of total cases) 
or 10-1-10 M0 total mass (0.2% of total cases). Orbital 
separations varied in the range 3-7 AU, with an average 
value of 5 AU (Table HI). 

In the case that two of the massive objects in the sys- 
tem fo rm a binary and source n is a high-mass YSO of 
8 Mq (iGreenhill et al. 2004b), the total mass of the clus- 
ter can be increased by swapping one of the 0.5 Mq stars 
for an 8 Mq star, keeping the same stellar density (total 
mass of 52.5 Mq). This would describe a system where 
Mi = 10 10 Mq, Mbn = 10 Mq, and M„ = 8 Mq. 
Simulations of the new system resulted in a similar num- 
ber of total runaway stars (19%) of maximum 8-10 Mq in 
0.7% of cases. Of these, one case out of 7 (0.1% of to- 
tal) gave runaways without formation of a binary. New 
binaries had total masses of 10-1-5 Mq in one case (0.1% 
of total), 10-1-8 Mq in one case (0.1% of total), and 
10+10 Mq in four cases (0.4% of total). Orbital sep- 
arations varied in the range 4-11 AU, with an average 
value of 6.5 AU. Since results are qualitatively similar to 
the system with total mass of 45 Mq, we did not report 
the outcome of the simulations for this case in Table |4l 

Our simulations show that, without the inclusion of 
very massive objects (^20 Mq) in the cluster, the prob- 
ability of ejecting massive members (> 8 Mq) from 
the system is very low (< 1 %). We note that pre- 
vious N-body simulations from iGomez et al.l (|2008f ) re- 
sulted in a larger fraction of ejections of objects with 
8 Mq (15%) because they assumed several cluster mem- 
bers with masses in the range 8-20 Mq and allowed for 



formation of binaries as massive as 36 Mq. However, 
such massive objects cannot readily be identified with 
any of the known objects in the BN/KL region. 

6.4. Caveat with the decay of a multiple system of single 
stars: Implications of binary formation on the disk 
around Source I 

As argued from the results of N-body simulations for 
systems of three or more stars. Source I is likely to be a 
binary comprising stars of comparably high mass. The 
formation of a tight binary implies a close passage (0(10 
AU) or less) , well within the uncertainties in the analysis 
of proper motion data (see Sect. [S]). However, the exis- 
tence of disks or other structures in close proximity to BN 
and Source I poses a challenge to explain. We observe 
an ionized disk with a ^ 50 AU radius around Source I, 
while BN has associated dust and (molecular and ion- 
ized) g as material in a disk-like structure of tens of AU in 
radius (| Jiang et al.l[2005t iRodnguez et al.l 120091 ). Close 
passages truncate disks larger in radius than the mini- 
mum stellar separation or destroy them all together. In 
the following, we discuss two possibilities: 1) the original 
disks are destroyed in the encounter and then rebuilt in 
the following 500 years; 2) the disks are truncated but 
survive the encounter. 

In the first hypothesis, evidence of a disk around 
Source I today is indicative of a viable rebuilding mecha- 
nism fueled by cither residual material dispersed from the 
original disks or interstellar gas in the region. We con- 
sider two possible mechanisms for rebuilding the disk: 
Bondi-Hoyle (B-H) accretion and tidal disruption of a 
low-mass object. 

BH accretion is the physical mechanism through which 
a compact object passing throu gh an ambient medium 
gains mass ()Bondi &: Hovlelll944i ) . Material within a crit- 
ical radius (BH radius) becomes bound and accreted by 
the protostar. Equating the potential and kinetic en- 
ergies, we define the BH radius Rbh = GM/v^, an 
associated accretion rate Mbh = t^pRbh''^' accre- 
tion timescale tBH ~ Rbh/v, and a disk mass Md = 
Mbh /tBH , where M is the compact-object mass and 

V is the velocity relative to the surrounding medium 
of density p. Assuming M = 20 Mq (Sect. 16. 3|) and 

V = 12 km s~^for Source I (Sect. E]), and n{H2) = 
10^ cm^'^ for the molecular gas in the region (as in the 
hot core), we obtain: Rbh = 120 AU, tBH = 66 years, 
Mbh = 5 X 10"^ Mq yr^^, M^ = 0.002 Mq. Although 
the radius estimate is consistent with the observations, 
a more sophisticated modeling of BH flow, using 2D and 
3D s imulations for both subsonic and turbulent flows 
(e.g. iKrumholz et al.l [2005t lEdgad [200I ) . demonstrates 
that gravitational focusing occurs from the initial scale 
of Rbh to the scale of an irregular accr etion pseudo-disk 
of radius < 0.1 x Rbh (e.g. see Fig. 1 of lThroop fc Ballvl 
I2OO8D . Hence, it appears unlikely that BH accretion 
would support a disk with comparable radius and sym- 
metry to that observed. In addition, iMatthews et al.l 
(|2010l ) report a lower-limit of ~ 0.002 Mq to the disk 
mass, as estim ated from SiO maser emission, whereas 
iBeuther et al.l ([2004) report an upper limit of 0.2 Mq es- 
timated from the continuum flux density at 345 GHz, 
indicating that the mass gathered in 500 years via BH 
would be quite low. We conclude that BH accretion can- 
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not explain the symmetric morphology, the size, and the 
mass of the disk we observe around Source I. 

Another possible mechanism for building the disk 

around Sour ce I is suggested by simulations from 

iDavies et alJ (j2006) , which show that a close encounter 
between 10 Mq and 3 Mq stars can lead to a tidal dis- 
ruption of the low-density 3 Mq star and subsequent 
formation of a massive disk around the massive star. 
If, in the encounter, one low-mass star were stripped 
away, this process would contribute to create a massive 
disk/envelope around Source I, with possible rapid accre- 
tion of a significant fraction of it onto the protostar. This 
scenario, however, faces two major challenges. First, in 
order for the tidal disruption to be effective, the two 
stars must pass each other to within roughly one stellar 
radius, and by implication, the stellar densi ty of the clus- 
ter m ust be quite high, 2 x 10^ Mq pc~^ (jPavies et al.l 
I2006D . Second, material following a disruption must be 
redistributed from a characteristic radius of a few tenths 
of an AU (orbital separation of the two stars) to the 
observed disk radius of tens of AU. The timescale for 
disk spreading by viscosity is t^, ^ r'^/i^, for disk ra- 
dius r and viscosity oc a (a is the viscosity parame- 
ter). We obtain viscous timescales ranging from 10^ yr 
{r = 10 A U, a = IQ-^) up to l O*^ yr (r = 100 AU, 
a = 10-3; iHollenbach et al.ll2000f ). Tidal disruption is 
hence an improbable mechanism for the reformation of 
the disk around Source I in 500 years. 

We conclude that the observation of a relatively large 
and symmetric disk around Source I is inconsistent with 
the formation of a binary and reformation of a circumbi- 
nary disk after 500 years. 

We consider then the alternative hypothesis where the 
original disks survive, although severely affected by the 
dynamical interaction (e.g., truncation, mass redistribu- 
tion, accretion burst, etc.). In order to assess conditions 
under which a circumstellar disk is retained in a dynam- 
ical encounter or close passage, one should include disk 
particles in the simulations, which is beyond the scope of 
this article. A follow-up paper (Moeckel et al. in prep.) 
will consider a dynamical interaction between a binary 
associated with a disk and a single star. The probability 
of disk survival in an interaction between a binary and 
a single is higher than a complex three-body interaction. 
The main result is that a significant fraction of the disk 
is retained in the case of either a clean and fast exchange 
of companion in the binary (one or two close passages) 
or in the case of survival of the binary, where the original 
circumbinary disk is truncated at the radius of the close 
passage with the single star. 

6.5. (Favored) Scenario 4- Three-body interaction: BN 
and a pre-existing Source I binary 

Ejection of massive stars and formation of tight bina- 
ries in multiple-body (>3) interactions have relatively 
low probability (< 1%) in the range of stellar den- 
sities and masses assumed in our N-body simulations 
(Sect. 16. 3p . Despite the fact that binaries are unlikely 
to arise by chance encounters, observations show that 
a large fraction (> 60%) of massive stars in OB clus- 
ters a nd/or associations are binary and multiple sys- 
tems ijMason et al.l 120091) . The abundance of bina- 
ries and triples inferred from observations suggests that 
most stars are born that way ("primordial" binaries). 



(a) 




Fig. 7. — Three body simulation for the interaction between 
a Source I binary (10+10 Mq; shown in red) and BN (10 Mq; 
shown in black) at two epochs. Initial positions for the simulations 
are marked by filled circles, with the Source I binary components 
separated by 10 AU and BN at a distance of 500 AU. Initial ve- 
locities are due to the thermal velocity of the region, vtjn,rm=0-4 
km s~^. Arrows mark source positions and 3D motion directions 
at times after closest Source I - BN approach of (a) tftiSO years 
and (b) tsd500 years. After 500 years, the Source I binary and BN 
have positions that are ~3200 AU apart and are moving away from 
each other with 3D velocities of 14 and 29 km s" ^ respectively, in 
agreement with the observations. 

Combined hydrodynamic and N-body simulations have 
shown that primordial binaries can form from both 
fragmentation of colla psing cores and accretion disks 
(| Alexander et al.ll2008l : iMoeckel fc Batell2010D . 

In the following, we reconsider a 3-body interaction 
where two stars lie in a (primordial) binary (Source I), 
and a close passage occurs by the third (BN). We as- 
sume a binary of total mass 20 Mq, with equal pri- 
mary and secondary masses of 10 Mq and an orbital 
separation of 10 AU. The third star has a mass of 
10 Mq and is at an initial distance of 500 AU. We 
have numerically integrated 1000 cases with three bod- 
ies in a bound system starting with a velocity dispersion 
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vtherm = 0.4 km s^^ and followed their dynamical inter- 
actions for 1000 years. Figure[7]shows the outcome of our 
simulations for one specific case. In 39% of cases, both 
the single and the binary stars acquire very high speeds 
in their interaction, flying apart in opposite directions 
with velocities in the range 10-90 km s~^ . In all cases, 
the single 10 Mq star runs away from the center of inter- 
action with approximately twice the velocity of the bi- 
nary, in agreement with linear momentum conservation. 
In the cases with ejections, the original binary survives 
in 15% of cases, while in 24% of cases a new binary is 
formed by exchanging one component (39% total). The 
orbital separation of binary components, averaged over 
the last 50 years of the simulation, varies in the range 
1-9 AU, with an average of 4 AU. The binding energy of 
the tighter binary (initially at a = 10 AU) compensates 
for the kinetic energy acquired by the two runaway stars, 
in agreement with conservation of mechanical energy. In 
the cases where the original binary survives (15%), the 
minimum periastron (over 1000 years) between the sin- 
gle and the binary varied in the range 1-60 AU, with an 
average value of 16 AU. 

The new 3-body simulations show that in a high frac- 
tion of cases, massive stars are ejected (39%) and the 
original binary survives (15%). Stellar ejections require a 
massive binary (~ 20 M0) as initial condition. In explor- 
ing parameter space, we also simulated a 10 Mq star in 
parabolic orbit about a binary of total mass 11 M©, with 
mass distributed in 10-f 1 M© or 8-1-3 M© units. With 
these starting conditions, only a small fraction (< 1%) 
resulted in ejection, and the star that was lost was the 
lower mass member. The process necessarily involved 
exchange of binary members, where the more massive 
companion and the third massive object formed an al- 
most equal-mass binary whereas the lower mass member 
was ejected from the system. At the moment, we favor 
the cases where the binary survives without exchanging 
components, because therein the inner reaches of the ac- 
cretion disk can be preserved (Sect. 16.41) . A follow-up 
paper (Moeckel et al. in prep.) will investigate the ef- 
fects of a stellar dynamical interaction on a circumbinary 
disk as a function of stellar mass, binary orbital separa- 
tion, single-binary periastron, and swapping of binary 
members. 

7. MASS OF SOURCE I 

The assumed dynamical interaction between Source I 
and BN, is concomitant wi th an estimated mass of ~ 
20 M© for Source I. Recently [Matthews et all ()2010[ ) esti- 
mated a dynamical mass using SiO maser proper motion 
data. This provided a lower limit of > 7 M©. The in- 
fluence of non-gravitational forces on gas dynamics may 
be resp onsible for the discrep ancy between dynamical 
masses. [Matthews et al.l ()2010D detected both radial and 
rotational motions in the inferred disk around the YSO, 
at radii < 100 AU. Although Doppler velocity data in 
the outflow show differential rotation of material, the 
actual rotation curve may be sub-Keplerian, e.g., due 
to radiative and/or magnetic forces, which can lead to 
underestimate s of the cent ral object mass. In particu- 
lar, lMatthcws ~et al.l ()2010f ) reported possibly bent tra- 
jectories among proper motions and suggest that mag- 
netic fields may play a role in . shapi ng the gas dynamics 
around Source I. lGirart et al.l ()2009l ) found recently that 



the evolution of the gravitational collapse of a massive 
hot molecular core in G31. 41-1-0. 31 is controlled by mag- 
netic fields (on scales ~ 10000 AU), which appear to 
be effective in removing angular momentum and slow- 
ing down circumst ellar gas rotation , as expecte d from 
magnetic braking (iGalli et al.l 120061 ). IShu et al.l (|2007f ) 
considered the case of an accretion disk threaded with 
magnetic fields that are squeezed in (scales < 1000 AU) 
from the interstellar medium by gravitational collapse. 
They find that the poloidal component of the magnetic 
field produces, as a consequence of its rotation, a change 
in the radial force balance in the disk, giving a deviation 
from a Keplerian profile. The resulting rotation law has 
the same power-law dependence with radius as a Keple- 
rian profile, but scaled with a coefficient / < 1, which in- 
cludes the effects of the magnetic field. As a consequence 
of this deviation from Keplerian rotation, the dynamical 
mass of the central object derived from an observed rota- 
tion curve will be smaller than its actual value of a factor 
Mobs/Mtrue = where Mtme is the true mass and Mobs 
is the observed mass. The departure from Keplerian ro- 
tation can be analytically expressed in terms of the mag- 
netic coefficient / as: I - p = 0.5444M,/(AgMdia3e) 
()Shu et al.l l2007f ) . where Md is the accretion rate onto 
the protostar with mass M,, tage is the stellar age, and 
Ao is a numerical coefficient from the model (assumed 
to be 4).^ Assuming for Source I Mtme ~ 20 MQ(from 
conservation of linear momentum of the BN-I system) 
and Mobs ^ 8 M© (from SiO maser dynamics) , we derive 
P ^ 0.4. Based on the formula above, this value can 
be readily achieved by assuming reasonable values for 
mass accretion rates and age: Md = 10^^ M© yr^^and 
tage = 10^ years (for Aq =4). This analysis shows that, 
in the assumption that magnetic fields support the disk, 
more than half of the mass of the central object can be 
easily hidden from observations. Nevertheless, we cau- 
tion that the estimated parameters are indicative and 
that careful modeling is required to determine the true 
mass of Source I. Future direct measurements of mag- 
netic field amplitude in the disk around Source I, and an 
accurate estimate of the mass-accretion rate as well as 
the stellar age, will enable a better understanding of the 
effects of magnetic fields on mass measurement from gas 
dynamics. 

8. WHAT IS POWERING THE WIDE-ANGLE H2 OUTFLOW 
IN BN/KL? 

The origin, nature, and powering source of the fast, 
poorly collimated bipolar outfiow in BN/KL have been 
long debated. The H2 finger structure consists of over 
100 individual bow shocks with a large opening angle 
and indicates tha t a powerful eveiit occu rred in the center 
of BN/KL (e.g., lAllen fc BurtOTllll993[) . Proper motion 
studies of Herbig-Haro (HH) objects in BN/KL from the 
ground (iLee & Burton"2000') and with the Hubble Space 
Telescope (Doi et al. 2002) showed that more distant 
bullets are moving faster and indicate an origin about 
1000 y ears ago (assuming no deceleration) . IZapata et al.l 
(j2009lf ) observed the CO (2-1) line with the Submillimeter 

* We note that this analysis assumes magnetic fields thread the 
disk and ignores the stellar magnetosphere. The formula above is 
then generally applicable to different kinds of stars, i.e. convective 
(magnetic) and radiative (non- magnetic), or, equivalently, low- and 
high-mass stars, respectively. 
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Array and identified several filaments arising from a com- 
mon center and showing linear velocit y grad ients with 
distance from the center. IZapata et alj ()2009t ) proposed 
that the outflow was produced by the disintegration of 
the young protostellar cluster formed by BN, I, and n. 

We propose a similar scenario, where the dynamical 
interaction between BN and the Source I binary, as de- 
scribed in Sect. 16. 5[ may have produced the fast CO out- 
flow and associated H2 fingers in BN/KL. The outfiow 
has a mass of about 10 Mc?, and a kinetic energy of about 
3 X lO''^ ergs (iKwan fc Scovillel[T976il . The hardening of 
a 20 M© binary to an orbital separation of 2 AU would 
release about 5 x 10'''^ ergs of gravitational potential en- 
ergy, enough to drive the powerful outfiow (3 x lO''^ ergs) 
and compensate for the kinetic energy of the BN-I system 
(10^^ erg). This outfiow would not be the typical outfiow 
powered by accretion onto a protostar, but rather the 
relic of a one-time event. In this picture, the expanding 
material was ejected at about the same time and acceler- 
ated from the area where the BN-I encounter took place, 
in agreement with the general "Hubble-flow" trend ev- 
idenced by proper motio ns of HH objects a nd velocity 
gradients in CO emission (jZapata et al.l[2009| ). Since the 
disrupted circumstellar material is dispersed in all direc- 
tions, this picture explains why the present structure of 
Source I is not related to the geometry of the outflow on 
large scales. 

We note that, although the formation or hardening 
of a massive binary in a dynamical interaction can ex- 
plain energetically the expansion of the outflowing gas, 
the physical mechanism at its origin remains unknown. 
More sophisticated combined N-body and hydrodynamic 
simulations that include also material in circumstellar 
disks and envelopes may help to understand the fate of 
circumstellar gas associated with individual protostars 
participating in a dynamical encounter. 

9. DISCUSSION AND CONCLUSIONS 

We have presented high-sensitivity, multi-epoch VLA 
A7 mm continuum observations of the radio sources in 
the Orion BN/KL star-forming region. These observa- 
tions have shown that the radio emission from Source I 
is elongated NW-SE and is stable over a decade. The 
morphology and the stability of the emission is consis- 
tent with an ionized edge-on disk. 

We measured the proper motions of Source I and the 
BN object for the flrst time at 43 GHz, confirming that 
they are moving with high speeds (12 and 26 km s^"'^ , 
respectively) approximately in opposite direction with 
respect to each other. We discussed possible dynamical 
scenarios that can explain the anomalously large motions 
of both BN and Source I, and implications for the dynam- 
ics of the whole BN/KL region. Our new measurements 
are inconsistent with a close encount er between BN an d 
Ori C in the past, as proposed bvlTanI (j2004 l2008b[ ). 
They confirm the previously proposed scenario of a dy- 
namical interactio n between Sou r ce I a nd BN 500 years 
ago. In particular, iGomez et al.l ()2008[ ) proposed that a 
dynamical decay of a protostellar cluster resulted in the 
formation of a tight binary (Source I) and the ejection of 
other two sources (BN and n). However, the presence of 
a relatively large, massive disk around Source I and the 
low-probability of a multiple-body interaction among sin- 
gle stars, excludes that Source I formed a binary in the 



course of the encounter and subsequently managed to 
rebuild its disk in 500 years. We propose that Source I 
was a binary system before the close passage with BN, 
which possibly truncated the original circumbinary disk 
to ~ 50 AU. Inference that this is the natal accretion disk 
is in agreement with dis k-mediated accretion fo r Source I, 
as recentlv proposed bv lMatthews et al.l ()2010[ ). This dy- 
namical event would have resulted in the ejection of both 
stars and the hardening of the original binary. 

The dynamical scenarios discussed in this paper were 
tested with N-body numerical simulations of dynami- 
cal interactions among members of a protostellar cluster. 
These simulations, however, illustrate a simplified phys- 
ical process of gravitational interaction among point- 
masses. A follow-up study will include the effects of a 
circumbinary disk in the dynamical interactions between 
protostars. 

The mass of Source I is still controversial. From linear 
momentum conservation of the BN-I system, we estimate 
^ 20 Mq — much higher than the value e stimated from 
SiO maser dynamics, > 7 M© CMatthews et al. 2010). On 
one hand, a strong argument in support of the dynami- 
cal estimate is that two independent families of motions, 
Keplerian rotation in a disk and escape speed in an out- 
flow, provide similar mass values, although relevant to 
different radii and force balance. On the other hand, the 
analysis in this paper is strongly based on the fact that 
5 to 6 observables (proper motions, LOS velocities, and 
positions in the plane of the sky) are consistent with the 
hypothesis of a close passage between Source I and BN. 
One possibility for reconciling the two models is that non- 
gravitational forces may play a significant role in driving 
the dynamics of molecular gas around Source I. Assum- 
ing the disk is magnetically supported, a significant frac- 
tion of the mass of the central object can be hidden in the 
Keplerian velocity profile of maser emission, leading to 
an underestimate of the mass (|Shu et al.ll20"07h . We note 
that an equal-mass binary with total mass 10 M© would 
not produce an HII region, and an unequal-mass bi- 
nary with the same total mass would have a very low- 
probability of surviving in a dynamical interaction with 
an equally-massive body (e.g.. Sect. [6?5|) . Hence, unless 
we assume the presence of another very massive source 
in the region (e.g., n), the dynamic scenario proposed 
here suggests a central mass ~ 20 Mq for Source I. We 
explicitly note that, owing to the non-linear dependence 
of luminosity on mass, the total luminosity developed by 
a binary composed of two 10 M© stars cannot explain 
the high luminosity {'^ 10^ Lq) of the BN/KL nebula. 
The gravitational energy released by the hardening of a 
binary with total mass 20 M© would, however, be suf- 
ficient to power the gas expansion observed in H2 and 
CO emission in the BN/KL nebula, although a detailed 
physical mechanism for such a process is still lacking. 

A final, general remark concerns models of massive star 
formation. Source I is probably the best case known of 
a ma ssive protostar w ith ongoing disk-mediated accre- 
tion (|Matthews et al.l 12010) . This work shows that re- 
cent dynamical interactions played a fundamental role 
in shaping the protostellar properties in Orion BN/KL. 
Star formation in this region might be then atypical in 
the context of the classic theory of isolated low-mass star 
formation (Shu et al.. .19871) . Nevertheless, considering 
that most (although not all) massive stars are formed in 
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dense protostellar clusters ()Lada fc Ladall20d3l ). dynam- 
ical interactions in very early stages of cluster evolution 
may be common in the context of crowded, massive star 
formation. Proper motion studies of the radio/mm con- 
tinuum sources in similar regions are required to confirm 
this hypothesis. The EVLA, when fully commissioned, 
will offer broadband (up to 8 GHz) imaging, resulting 
in a factor of 10 improvement in continuum sensitivity, 
making this kind of study possible in principle for regions 
three times more distant. ALMA, with a similar or even 
better resolution, will enable similar studies at (sub)mm 
wavelengths. Systematic studies with these instruments 



will allow us to assess the role of dynamical interactions 
in massive star formation. 
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APPENDIX 
ABSOLUTE ASTROMETRIC ACCURACY 

The absolute astrometry in each of the program datasets reported here has been determined with high accuracy 
(of order ^^5 mas). The error at each epoch is taken to be the square root of the sum in quadrature of the 4 main 
contributions to the total error budget: 

Cp = ^e1 + el + el + el 

where individual contributions account for angular separation (calibrator to target), noise limited or fitting uncertainty, 
source structure, and frequency-dependent errors. Below we describe in detail individual error contributions: 

• et is the theoretical error in computing the absolute astrometry from a calibrator with a separation of 1° 
from the target. For a typical baseline length of 10 km and accuracy of 1 cm (appropriate for the VLA in 
A configuration), and a calibrator-target separation of 1.4°, one derives et ~ 5 mas. This is the dominant 
contribution to the overall error budget. 

• e„ is the theoretical, noise- limited positional uncertainty, given by e„ = 0.5 9/SNR, where 9 is the FWHM of 
the synthesized beamwidth of the array, and SNR is the peak intensity divided by the RMS noise in the map. 
We compared this error with the fitting error from JMFIT and took the larger of the two (typically ~ 1 mas). 

• is the uncertainty introduced by the source structure. Since both BN and Source I are resolved at 7 mm 
in A-configuration, different weighting schemes result in slight changes in morphology, and possibly in different 
peak positions. We quantified this contribution by taking the dispersion (i.e., max separation) of the positions 
in the images with different weighting and restoring beams: ~1 mas for BN and I (resolved sources) and 
Bs ~0 mas for H and n (point sources). 

• is the error introduced by the application of the self-calibration solutions from the maser to the continuum 
data, having a frequency separation Av ^ 100 — 350 MHz. In previous VLA experiments, we es tablished an 
angular offset per MHz as large as 0.01 mas/MHz introduced by calibrating one band by another (|Goddi et al.l 
[2009). We cannot quantify this contribution for the programs described here (we observed only at one frequency 
offset from the v = 1 maser line), but assuming this contribution does not change much from program to program, 
one can anticipate a contribution at most of a few mas. 

We report typical values of individual contributions in Table 2. By adding up all contributions, typical total 
uncertainties are ^ 5 mas. 
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